﻿
namespace GMap.NET
{
    using System;
    using System.Collections.Generic;
    using System.Diagnostics;

    /// <summary>
    /// defines projection
    /// </summary>
    public abstract class PureProjection
    {
        readonly List<Dictionary<PointLatLng, GPoint>> FromLatLngToPixelCache = new List<Dictionary<PointLatLng, GPoint>>(33);
        readonly List<Dictionary<GPoint, PointLatLng>> FromPixelToLatLngCache = new List<Dictionary<GPoint, PointLatLng>>(33);

        public PureProjection()
        {
            for (int i = 0; i < FromLatLngToPixelCache.Capacity; i++)
            {
                FromLatLngToPixelCache.Add(new Dictionary<PointLatLng, GPoint>());
                FromPixelToLatLngCache.Add(new Dictionary<GPoint, PointLatLng>());
            }
        }

        /// <summary>
        /// size of tile
        /// </summary>
        public abstract GSize TileSize
        {
            get;
        }

        /// <summary>
        /// Semi-major axis of ellipsoid, in meters
        /// </summary>
        public abstract double Axis
        {
            get;
        }

        /// <summary>
        /// Flattening of ellipsoid
        /// </summary>
        public abstract double Flattening
        {
            get;
        }

        /// <summary>
        /// get pixel coordinates from lat/lng
        /// </summary>
        /// <param name="lat"></param>
        /// <param name="lng"></param>
        /// <param name="zoom"></param>
        /// <returns></returns>
        public abstract GPoint FromLatLngToPixel(double lat, double lng, int zoom);

        /// <summary>
        /// gets lat/lng coordinates from pixel coordinates
        /// </summary>
        /// <param name="x"></param>
        /// <param name="y"></param>
        /// <param name="zoom"></param>
        /// <returns></returns>
        public abstract PointLatLng FromPixelToLatLng(long x, long y, int zoom);

        public GPoint FromLatLngToPixel(PointLatLng p, int zoom)
        {
            return FromLatLngToPixel(p, zoom, false);
        }

        /// <summary>
        /// get pixel coordinates from lat/lng
        /// </summary>
        /// <param name="p"></param>
        /// <param name="zoom"></param>
        /// <returns></returns>
        public GPoint FromLatLngToPixel(PointLatLng p, int zoom, bool useCache)
        {
            if (useCache)
            {
                GPoint ret = GPoint.Empty;
                if (!FromLatLngToPixelCache[zoom].TryGetValue(p, out ret))
                {
                    ret = FromLatLngToPixel(p.Lat, p.Lng, zoom);
                    FromLatLngToPixelCache[zoom].Add(p, ret);

                    // for reverse cache
                    if (!FromPixelToLatLngCache[zoom].ContainsKey(ret))
                    {
                        FromPixelToLatLngCache[zoom].Add(ret, p);
                    }

                    Debug.WriteLine("FromLatLngToPixelCache[" + zoom + "] added " + p + " with " + ret);
                }
                return ret;
            }
            else
            {
                return FromLatLngToPixel(p.Lat, p.Lng, zoom);
            }
        }

        public PointLatLng FromPixelToLatLng(GPoint p, int zoom)
        {
            return FromPixelToLatLng(p, zoom, false);
        }

        /// <summary>
        /// gets lat/lng coordinates from pixel coordinates
        /// </summary>
        /// <param name="p"></param>
        /// <param name="zoom"></param>
        /// <returns></returns>
        public PointLatLng FromPixelToLatLng(GPoint p, int zoom, bool useCache)
        {
            if (useCache)
            {
                PointLatLng ret = PointLatLng.Empty;
                if (!FromPixelToLatLngCache[zoom].TryGetValue(p, out ret))
                {
                    ret = FromPixelToLatLng(p.X, p.Y, zoom);
                    FromPixelToLatLngCache[zoom].Add(p, ret);

                    // for reverse cache
                    if (!FromLatLngToPixelCache[zoom].ContainsKey(ret))
                    {
                        FromLatLngToPixelCache[zoom].Add(ret, p);
                    }

                    Debug.WriteLine("FromPixelToLatLngCache[" + zoom + "] added " + p + " with " + ret);
                }
                return ret;
            }
            else
            {
                return FromPixelToLatLng(p.X, p.Y, zoom);
            }
        }

        /// <summary>
        /// gets tile coorddinate from pixel coordinates
        /// </summary>
        /// <param name="p"></param>
        /// <returns></returns>
        public virtual GPoint FromPixelToTileXY(GPoint p)
        {
            return new GPoint((long)(p.X / TileSize.Width), (long)(p.Y / TileSize.Height));
        }

        /// <summary>
        /// gets pixel coordinate from tile coordinate
        /// </summary>
        /// <param name="p"></param>
        /// <returns></returns>
        public virtual GPoint FromTileXYToPixel(GPoint p)
        {
            return new GPoint((p.X * TileSize.Width), (p.Y * TileSize.Height));
        }

        /// <summary>
        /// min. tile in tiles at custom zoom level
        /// </summary>
        /// <param name="zoom"></param>
        /// <returns></returns>
        public abstract GSize GetTileMatrixMinXY(int zoom);

        /// <summary>
        /// max. tile in tiles at custom zoom level
        /// </summary>
        /// <param name="zoom"></param>
        /// <returns></returns>
        public abstract GSize GetTileMatrixMaxXY(int zoom);

        /// <summary>
        /// gets matrix size in tiles
        /// </summary>
        /// <param name="zoom"></param>
        /// <returns></returns>
        public virtual GSize GetTileMatrixSizeXY(int zoom)
        {
            GSize sMin = GetTileMatrixMinXY(zoom);
            GSize sMax = GetTileMatrixMaxXY(zoom);

            return new GSize(sMax.Width - sMin.Width + 1, sMax.Height - sMin.Height + 1);
        }

        /// <summary>
        /// tile matrix size in pixels at custom zoom level
        /// </summary>
        /// <param name="zoom"></param>
        /// <returns></returns>
        public long GetTileMatrixItemCount(int zoom)
        {
            GSize s = GetTileMatrixSizeXY(zoom);
            return (s.Width * s.Height);
        }

        /// <summary>
        /// gets matrix size in pixels
        /// </summary>
        /// <param name="zoom"></param>
        /// <returns></returns>
        public virtual GSize GetTileMatrixSizePixel(int zoom)
        {
            GSize s = GetTileMatrixSizeXY(zoom);
            return new GSize(s.Width * TileSize.Width, s.Height * TileSize.Height);
        }

        //private static readonly object syncRoot = new object();

        /// <summary>
        /// gets all tiles in rect at specific zoom
        /// </summary>
        public List<GPoint> GetAreaTileList(RectLatLng rect, int zoom, int padding)
        {
            List<GPoint> ret = new List<GPoint>();

            GPoint topLeft = FromPixelToTileXY(FromLatLngToPixel(rect.LocationTopLeft, zoom));
            GPoint rightBottom = FromPixelToTileXY(FromLatLngToPixel(rect.LocationRightBottom, zoom));

            long begin_x = topLeft.X - padding;
            long end_x = rightBottom.X + padding;
            long begin_y = topLeft.Y - padding;
            long end_y = rightBottom.Y + padding;

            for (long x = begin_x; x <= end_x; ++x)
            {
                for (long y = begin_y; y <= end_y; ++y)
                {
                    GPoint p = new GPoint(x, y);
                    if (!ret.Contains(p) && p.X >= 0 && p.Y >= 0)
                    {
                        ret.Add(p);
                    }
                }
            }

            return ret;
        }

        public List<GPoint> GetAreaTileListBetweenZoom(RectLatLng rect, int minZoom, int maxZoom, int padding)
        {
            List<GPoint> ret = new List<GPoint>();
            for (int i = minZoom; i <= maxZoom; ++i)
            {
                ret.AddRange(GetAreaTileList(rect, i, padding));
            }

            return ret;
        }

        /// <summary>
        /// The ground resolution indicates the distance (in meters) on the ground that’s represented by a single pixel in the map.
        /// For example, at a ground resolution of 10 meters/pixel, each pixel represents a ground distance of 10 meters.
        /// </summary>
        /// <param name="zoom"></param>
        /// <param name="latitude"></param>
        /// <returns></returns>
        public virtual double GetGroundResolution(int zoom, double latitude)
        {
            return (Math.Cos(latitude * (Math.PI / 180)) * 2 * Math.PI * Axis) / GetTileMatrixSizePixel(zoom).Width;
        }

        /// <summary>
        /// gets boundaries
        /// </summary>
        public virtual RectLatLng Bounds
        {
            get
            {
                return RectLatLng.FromLTRB(-180, 90, 180, -90);
            }
        }

        public virtual PointLatLng TileOrigin
        {
            get
            {
                return new PointLatLng(this.Bounds.Lat, this.Bounds.Lng);
            }
        }

        public GPoint FromLatLngToTileXY(PointLatLng pointLatLng, int zoom)
        {
            return this.FromPixelToTileXY(this.FromLatLngToPixel(pointLatLng, zoom));
        }

        public GPoint FromLatLngToTileXY(double lat, double lng, int zoom)
        {
            GPoint p = this.FromLatLngToPixel(lat, lng, zoom);
            return this.FromPixelToTileXY(p);
        }

        public double DPI = 96.0;

        public int EpsgCode = 3857;

        /// <summary>
        /// Get level resolution
        /// </summary>
        /// <param name="level"></param>
        /// <returns></returns>
        public virtual double GetLevelResolution(int level)
        {
            return (((3.1415926535897931 * this.Axis) * 2.0) / (Math.Pow(2.0, (double)level) * this.TileSize.Width));
        }

        /// <summary>
        /// Get level scale
        /// </summary>
        /// <param name="level"></param>
        /// <returns></returns>
        public virtual double GetLevelScale(int level)
        {
            //if (this.Unit == MapUnit.Meter)
            if (this.EpsgCode == 3857)
            {
                return ((this.DPI / 0.0254) * this.GetLevelResolution(level));
            }
            return (((((this.GetLevelResolution(level) * this.DPI) * 2.0) * 3.1415926535897931) * this.Axis) / 9.144);
        }

        //public virtual double GetLevelScale(int level, double dpi)
        //{
        //    //if (this.Unit == MapUnit.Meter)
        //    if (this.EpsgCode == 3857)
        //    {
        //        return ((dpi / 0.0254) * this.GetLevelResolution(level));
        //    }
        //    return (((((this.GetLevelResolution(level) * dpi) * 2.0) * 3.1415926535897931) * this.Axis) / 9.144);
        //}

        public PointLatLng GetProjectedPoint(PointLatLng pointLatLng)
        {
            if (this.EpsgCode == 3857)
            {
                return LonLatToMercator(pointLatLng.Lng, pointLatLng.Lat);
            }
            //return pointLatLng.TransformFromWGS84(ProjectionUtil.GetWKTFromEpsgCode(this.EpsgCode));
            //EpsgCode=4326
            return new PointLatLng(pointLatLng.Lat, pointLatLng.Lng);
        }

        public PointLatLng LonLatToMercator(double x, double y)
        {
            double lng = ((x * 3.1415926535897931) * 6378137.0) / 180.0;
            double num = Math.Log(Math.Tan(((90.0 + y) * 3.1415926535897931) / 360.0)) / 0.017453292519943295;
            double lat = ((num * 3.1415926535897931) * 6378137.0) / 180.0;
            return new PointLatLng(lat, lng);
        }

        public PointLatLng MercatorToLonLat(double x, double y)
        {
            double lng = (x / (3.1415926535897931 * this.Axis)) * 180.0;
            y = (y / (3.1415926535897931 * this.Axis)) * 180.0;
            return new PointLatLng(57.295779513082323 * ((2.0 * Math.Atan(Math.Exp((y * 3.1415926535897931) / 180.0))) - 1.5707963267948966), lng);
        }

        #region -- math functions --

        /// <summary>
        /// PI
        /// </summary>
        protected static readonly double PI = Math.PI;

        /// <summary>
        /// Half of PI
        /// </summary>
        protected static readonly double HALF_PI = (PI * 0.5);

        /// <summary>
        /// PI * 2
        /// </summary>
        protected static readonly double TWO_PI = (PI * 2.0);

        /// <summary>
        /// EPSLoN
        /// </summary>
        protected static readonly double EPSLoN = 1.0e-10;

        /// <summary>
        /// MAX_VAL
        /// </summary>
        protected const double MAX_VAL = 4;

        /// <summary>
        /// MAXLONG
        /// </summary>
        protected static readonly double MAXLONG = 2147483647;

        /// <summary>
        /// DBLLONG
        /// </summary>
        protected static readonly double DBLLONG = 4.61168601e18;

        static readonly double R2D = 180 / Math.PI;
        static readonly double D2R = Math.PI / 180;

        public static double DegreesToRadians(double deg)
        {
            return (D2R * deg);
        }

        public static double RadiansToDegrees(double rad)
        {
            return (R2D * rad);
        }

        ///<summary>
        /// return the sign of an argument 
        ///</summary>
        protected static double Sign(double x)
        {
            if (x < 0.0)
                return (-1);
            else
                return (1);
        }

        protected static double AdjustLongitude(double x)
        {
            long count = 0;
            while (true)
            {
                if (Math.Abs(x) <= PI)
                    break;
                else
                   if (((long)Math.Abs(x / Math.PI)) < 2)
                    x = x - (Sign(x) * TWO_PI);
                else
                      if (((long)Math.Abs(x / TWO_PI)) < MAXLONG)
                {
                    x = x - (((long)(x / TWO_PI)) * TWO_PI);
                }
                else
                         if (((long)Math.Abs(x / (MAXLONG * TWO_PI))) < MAXLONG)
                {
                    x = x - (((long)(x / (MAXLONG * TWO_PI))) * (TWO_PI * MAXLONG));
                }
                else
                            if (((long)Math.Abs(x / (DBLLONG * TWO_PI))) < MAXLONG)
                {
                    x = x - (((long)(x / (DBLLONG * TWO_PI))) * (TWO_PI * DBLLONG));
                }
                else
                    x = x - (Sign(x) * TWO_PI);
                count++;
                if (count > MAX_VAL)
                    break;
            }
            return (x);
        }

        /// <summary>
        /// calculates the sine and cosine
        /// </summary>
        protected static void SinCos(double val, out double sin, out double cos)
        {
            sin = Math.Sin(val);
            cos = Math.Cos(val);
        }

        /// <summary>
        /// computes the constants e0, e1, e2, and e3 which are used
        /// in a series for calculating the distance along a meridian.
        /// </summary>
        /// <param name="x">represents the eccentricity squared</param>
        /// <returns></returns>
        protected static double e0fn(double x)
        {
            return (1.0 - 0.25 * x * (1.0 + x / 16.0 * (3.0 + 1.25 * x)));
        }

        protected static double e1fn(double x)
        {
            return (0.375 * x * (1.0 + 0.25 * x * (1.0 + 0.46875 * x)));
        }

        protected static double e2fn(double x)
        {
            return (0.05859375 * x * x * (1.0 + 0.75 * x));
        }

        protected static double e3fn(double x)
        {
            return (x * x * x * (35.0 / 3072.0));
        }

        /// <summary>
        /// computes the value of M which is the distance along a meridian
        /// from the Equator to latitude phi.
        /// </summary>
        protected static double mlfn(double e0, double e1, double e2, double e3, double phi)
        {
            return (e0 * phi - e1 * Math.Sin(2.0 * phi) + e2 * Math.Sin(4.0 * phi) - e3 * Math.Sin(6.0 * phi));
        }

        /// <summary>
        /// calculates UTM zone number
        /// </summary>
        /// <param name="lon">Longitude in degrees</param>
        /// <returns></returns>
        protected static long GetUTMzone(double lon)
        {
            return ((long)(((lon + 180.0) / 6.0) + 1.0));
        }

        /// <summary>
        /// Clips a number to the specified minimum and maximum values.
        /// </summary>
        /// <param name="n">The number to clip.</param>
        /// <param name="minValue">Minimum allowable value.</param>
        /// <param name="maxValue">Maximum allowable value.</param>
        /// <returns>The clipped value.</returns>
        protected static double Clip(double n, double minValue, double maxValue)
        {
            return Math.Min(Math.Max(n, minValue), maxValue);
        }

        /// <summary>
        /// distance (in km) between two points specified by latitude/longitude
        /// The Haversine formula, http://www.movable-type.co.uk/scripts/latlong.html
        /// </summary>
        /// <param name="p1"></param>
        /// <param name="p2"></param>
        /// <returns></returns>
        public double GetDistance(PointLatLng p1, PointLatLng p2)
        {
            double dLat1InRad = p1.Lat * (Math.PI / 180);
            double dLong1InRad = p1.Lng * (Math.PI / 180);
            double dLat2InRad = p2.Lat * (Math.PI / 180);
            double dLong2InRad = p2.Lng * (Math.PI / 180);
            double dLongitude = dLong2InRad - dLong1InRad;
            double dLatitude = dLat2InRad - dLat1InRad;
            double a = Math.Pow(Math.Sin(dLatitude / 2), 2) + Math.Cos(dLat1InRad) * Math.Cos(dLat2InRad) * Math.Pow(Math.Sin(dLongitude / 2), 2);
            double c = 2 * Math.Atan2(Math.Sqrt(a), Math.Sqrt(1 - a));
            double dDistance = (Axis / 1000.0) * c;
            return dDistance;
        }

        public double GetDistanceInPixels(GPoint point1, GPoint point2)
        {
            double a = (double)(point2.X - point1.X);
            double b = (double)(point2.Y - point1.Y);

            return Math.Sqrt(a * a + b * b);
        }

        /// <summary>
        /// Accepts two coordinates in degrees.
        /// </summary>
        /// <returns>A double value in degrees. From 0 to 360.</returns>
        public double GetBearing(PointLatLng p1, PointLatLng p2)
        {
            var latitude1 = DegreesToRadians(p1.Lat);
            var latitude2 = DegreesToRadians(p2.Lat);
            var longitudeDifference = DegreesToRadians(p2.Lng - p1.Lng);

            var y = Math.Sin(longitudeDifference) * Math.Cos(latitude2);
            var x = Math.Cos(latitude1) * Math.Sin(latitude2) - Math.Sin(latitude1) * Math.Cos(latitude2) * Math.Cos(longitudeDifference);

            return (RadiansToDegrees(Math.Atan2(y, x)) + 360) % 360;
        }

        /// <summary>
        /// Conversion from cartesian earth-sentered coordinates to geodetic coordinates in the given datum
        /// </summary>
        /// <param name="Lat"></param>
        /// <param name="Lon"></param>
        /// <param name="Height">Height above ellipsoid [m]</param>
        /// <param name="X"></param>
        /// <param name="Y"></param>
        /// <param name="Z"></param>
        public void FromGeodeticToCartesian(double Lat, double Lng, double Height, out double X, out double Y, out double Z)
        {
            Lat = (Math.PI / 180) * Lat;
            Lng = (Math.PI / 180) * Lng;

            double B = Axis * (1.0 - Flattening);
            double ee = 1.0 - (B / Axis) * (B / Axis);
            double N = (Axis / Math.Sqrt(1.0 - ee * Math.Sin(Lat) * Math.Sin(Lat)));

            X = (N + Height) * Math.Cos(Lat) * Math.Cos(Lng);
            Y = (N + Height) * Math.Cos(Lat) * Math.Sin(Lng);
            Z = (N * (B / Axis) * (B / Axis) + Height) * Math.Sin(Lat);
        }

        /// <summary>
        /// Conversion from cartesian earth-sentered coordinates to geodetic coordinates in the given datum
        /// </summary>
        /// <param name="X"></param>
        /// <param name="Y"></param>
        /// <param name="Z"></param>
        /// <param name="Lat"></param>
        /// <param name="Lon"></param>
        public void FromCartesianTGeodetic(double X, double Y, double Z, out double Lat, out double Lng)
        {
            double E = Flattening * (2.0 - Flattening);
            Lng = Math.Atan2(Y, X);

            double P = Math.Sqrt(X * X + Y * Y);
            double Theta = Math.Atan2(Z, (P * (1.0 - Flattening)));
            double st = Math.Sin(Theta);
            double ct = Math.Cos(Theta);
            Lat = Math.Atan2(Z + E / (1.0 - Flattening) * Axis * st * st * st, P - E * Axis * ct * ct * ct);

            Lat /= (Math.PI / 180);
            Lng /= (Math.PI / 180);
        }

        #endregion
    }
}
